Análise exploratória da base REDD, aplicando Classificação de Séries Temporais baseada em Aprendizado Supervisionado.
A base REDD contempla dados de consumo energético de 6 casas distintas. Primeiro, após a análise exploratória, será treinado e testado um modelo a partir dos dados da residência 1 e avaliado o quão bem esse modelo generaliza para os padrões não-observados da residência 2.
In [11]:
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
%matplotlib inline
from IPython.display import display
import datetime
import time, os
import math
import warnings
warnings.filterwarnings("ignore")
import glob
import nilmtk
In [2]:
PATH_DATASET = './datasets/REDD/'
def read_label():
label = {}
for i in range(1, 7):
hi = os.path.join(PATH_DATASET, 'low_freq/house_{}/labels.dat').format(i)
label[i] = {}
with open(hi) as f:
for line in f:
splitted_line = line.split(' ')
label[i][int(splitted_line[0])] = splitted_line[1].strip() + '_' + splitted_line[0]
return label
labels = read_label()
for i in range(1,3):
print('Residência {}: '.format(i), labels[i], '\n')
In [3]:
def read_merge_data(house):
path = os.path.join(PATH_DATASET, 'low_freq/house_{}/').format(house)
file = path + 'channel_1.dat'
df = pd.read_table(file, sep = ' ', names = ['unix_time', labels[house][1]],
dtype = {'unix_time': 'int64', labels[house][1]:'float64'})
num_apps = len(glob.glob(path + 'channel*'))
for i in range(2, num_apps + 1):
file = path + 'channel_{}.dat'.format(i)
data = pd.read_table(file, sep = ' ', names = ['unix_time', labels[house][i]],
dtype = {'unix_time': 'int64', labels[house][i]:'float64'})
df = pd.merge(df, data, how = 'inner', on = 'unix_time')
df['timestamp'] = df['unix_time'].astype("datetime64[s]")
df = df.set_index(df['timestamp'].values)
df.drop(['unix_time','timestamp'], axis=1, inplace=True)
return df
df = {}
for i in range(1,3):
df[i] = read_merge_data(i)
In [4]:
for i in range(1,3):
print('Shape dos dados da Residência {}: '.format(i), df[i].shape)
display(df[i].tail(3))
In [5]:
dates = {}
for i in range(1,3):
dates[i] = [str(time)[:10] for time in df[i].index.values]
dates[i] = sorted(list(set(dates[i])))
print('Os dados da Residência {0} contém medições de {1} dia(s) (de {2} a {3}).'.format(i,len(dates[i]),dates[i][0], dates[i][-1]))
print(dates[i], '\n')
In [6]:
# Plotar 2 primeiros dias de dados das residências 1 e 2
def plot_df(df, title):
apps = df.columns.values
num_apps = len(apps)
fig, axes = plt.subplots((num_apps+1)//2,2, figsize=(24, num_apps*2) )
for i, key in enumerate(apps):
axes.flat[i].plot(df[key], alpha = 0.6)
axes.flat[i].set_title(key, fontsize = '15')
plt.suptitle(title, fontsize = '30')
fig.tight_layout()
fig.subplots_adjust(top=0.95)
for i in range(1,3):
plot_df(df[i].ix[:dates[i][1]], 'Registros dos 2 primeiros dias da Residência {}'.format(i))
In [7]:
# Plot total energy sonsumption of each appliance from two houses
fig, axes = plt.subplots(1,2,figsize=(24, 10))
plt.suptitle('Energia total consumida por cada aparelho', fontsize = 30)
cons1 = df[1][df[1].columns.values[2:]].sum().sort_values(ascending=False)
app1 = cons1.index
y_pos1 = np.arange(len(app1))
axes[0].bar(y_pos1, cons1.values, alpha=0.6)
plt.sca(axes[0])
plt.xticks(y_pos1, app1, rotation = 90, fontsize=16)
plt.title('Residência 1')
cons2 = df[2][df[2].columns.values[2:]].sum().sort_values(ascending=False)
app2 = cons2.index
y_pos2 = np.arange(len(app2))
axes[1].bar(y_pos2, cons2.values, alpha=0.6)
plt.sca(axes[1])
plt.xticks(y_pos2, app2, rotation = 90, fontsize=16)
plt.title('Residência 2')
Out[7]:
In [8]:
# Split de treino, teste e validação
df1_train = df[1].ix[:dates[1][10]]
df1_val = df[1].ix[dates[1][11]:dates[1][16]]
df1_test = df[1].ix[dates[1][17]:]
print('df_train.shape: ', df1_train.shape)
print('df_val.shape: ', df1_val.shape)
print('df_test.shape: ', df1_test.shape)
In [9]:
# Exemplo de base, com X = ('mains_1','mains_2') e Y = (refrigerator_5)
# A previsão (desagregação) da corrente é na vertical
df_sample = df1_val[['mains_1','mains_2','refrigerator_5']]
df_sample.head(10)
Out[9]:
In [10]:
print('Dias compreendidos na leitura/desagregação:')
set([str(dt).split(' ')[0] for dt in df_sample.index])
Out[10]:
In [11]:
# Usando a corrente 1 e 2 (variaveis independetes) para a previsão do refrigerador (variavel dependente)
X_train1 = df1_train[['mains_1','mains_2']].values
y_train1 = df1_train['refrigerator_5'].values
X_val1 = df1_val[['mains_1','mains_2']].values
y_val1 = df1_val['refrigerator_5'].values
X_test1 = df1_test[['mains_1','mains_2']].values
y_test1 = df1_test['refrigerator_5'].values
print(
X_train1.shape, y_train1.shape,
X_val1.shape, y_val1.shape,
X_test1.shape, y_test1.shape
)
In [12]:
# Metrcas de avaliação da regressão
def mse_loss(y_predict, y):
return np.mean(np.square(y_predict - y))
def mae_loss(y_predict, y):
return np.mean(np.abs(y_predict - y))
# Serão usados os dados de validação para ajustar o parâmetros min_samples_split
min_samples_split = np.arange(2, 400, 10)
# Treinando o modelo
from sklearn.tree import DecisionTreeRegressor
def tree_reg(X_train, y_train, X_val, y_val, min_samples_split):
clfs = []
losses = []
start = time.time()
for split in min_samples_split:
clf = DecisionTreeRegressor(min_samples_split = split)
clf.fit(X_train, y_train)
y_predict_val = clf.predict(X_val)
clfs.append(clf)
losses.append( mse_loss(y_predict_val, y_val) )
print('Tempo de execução (s): ', round(time.time() - start, 0))
return clfs, losses
tree_clfs_1, tree_losses_1 = tree_reg(X_train1, y_train1, X_val1, y_val1, min_samples_split)
In [13]:
def plot_losses(losses, min_samples_split):
index = np.arange(len(min_samples_split))
bar_width = 0.4
opacity = 0.35
plt.bar(index, losses, bar_width, alpha=opacity, color='b')
plt.xlabel('min_samples_split', fontsize=30)
plt.ylabel('loss', fontsize=30)
plt.title('Loss (Validação) x min_samples_split', fontsize = '25')
plt.xticks(index + bar_width/2, min_samples_split, fontsize=20 )
plt.yticks(fontsize=20 )
plt.rcParams["figure.figsize"] = [24,15]
plt.tight_layout()
plot_losses(tree_losses_1, min_samples_split)
In [14]:
# Escolhendo o melhor modelo (minsplit x loss) e prevendo o consumo do refrigerador no conjunto de teste
ind = np.argmin(tree_losses_1)
tree_clf_1 = tree_clfs_1[ind]
y_test_predict_1 = tree_clf_1.predict(X_test1)
mse_tree_1 = mse_loss(y_test_predict_1, y_test1)
mae_tree_1 = mae_loss(y_test_predict_1, y_test1)
print('MSE no Conjunto de Teste:', mse_tree_1)
print('MAE no Conjunto de Teste:', mae_tree_1)
In [15]:
# Plotando os cnsumos REAL e o PREVISTO do refrigerador nos 6 dias dos dados de teste
def plot_each_app(df, dates, predict, y_test, title, look_back = 0):
num_date = len(dates)
fig, axes = plt.subplots(num_date,1,figsize=(24, num_date*5) )
plt.suptitle(title, fontsize = '25')
fig.tight_layout()
fig.subplots_adjust(top=0.95)
for i in range(num_date):
if i == 0: l = 0
ind = df.ix[dates[i]].index[look_back:]
axes.flat[i].plot(ind, y_test[l:l+len(ind)], color = 'blue', alpha = 0.6, label = 'REAL')
axes.flat[i].plot(ind, predict[l:l+len(ind)], color = 'red', alpha = 0.6, label = 'PREVISTO')
axes.flat[i].legend()
l = len(ind)
plot_each_app(df1_test, dates[1][17:], y_test_predict_1, y_test1, 'Consumo Real/Previsto do refrigerador nos 6 dias da Residência 1')
In [29]:
X_2 = df[2][['mains_2','mains_1']].values # Mesmas variáveis independentes, correntes 1 e 2)
y_2 = df[2]['refrigerator_9'].values
print(X_2.shape, y_2.shape)
In [30]:
y_predict_2 = tree_clf_1.predict(X_2)
mse_tree_2 = mse_loss(y_predict_2, y_2)
mae_tree_2 = mae_loss(y_predict_2, y_2)
print('MSE no Conjunto de Teste:', mse_tree_2)
print('MAE no Conjunto de Teste:', mae_tree_2)
In [31]:
plot_each_app(df[2], dates[2], y_predict_2, y_2, 'Modelo de Árvore de Decisão aplicado ao Refrigerador: treinado na Res. 1, prevendo na Res. 2')
In [36]:
# Lista de outros aparelhos da res. 1
appliances = list(df[1].columns.values[2:])
appliances.pop(2)
print(appliances)
In [37]:
# Treinando o modelo de Árvore para os outros aparelhos
def tree_reg_mult_apps():
start = time.time()
min_samples_split=np.arange(2, 400, 10)
pred = {}
for app in appliances:
list_clfs = []
losses = []
y_train = df1_train[app].values
y_val = df1_val[app].values
for split in min_samples_split:
clf = DecisionTreeRegressor(min_samples_split = split)
clf.fit(X_train1, y_train)
y_predict_val = clf.predict(X_val1)
list_clfs.append(clf)
losses.append( mse_loss(y_predict_val, y_val) )
ind = np.argmin(losses)
pred[app] = list_clfs[ind].predict(X_test1)
print('Tempo de execução (s): ', round(time.time() - start, 0))
return pred
mul_pred = tree_reg_mult_apps()
In [41]:
# Calculando erros/perdas (losses) de múltiplos aparelhos
def error_mul_app(mul_pred):
mse_losses = {}
mae_losses = {}
for app in appliances:
mse_losses[app] = mse_loss(mul_pred[app], df1_test[app].values)
mae_losses[app] = mae_loss(mul_pred[app], df1_test[app].values)
return mse_losses, mae_losses
mul_mse_tree, mul_mae_tree = error_mul_app(mul_pred)
In [42]:
for app in appliances:
m = np.mean(df1_test[app].values)
print('Consumo Médio de {0}: {1:.2f} - MSE: {2:.2f} - MAE: {3:.2f}'.format(app, m ,mul_mse_tree[app], mul_mae_tree[app]))
In [54]:
for app in appliances:
plot_each_app(df1_test, dates[1][17:], mul_pred[app], df1_test[app].values,
'{} - Real e Previsto nos 6 dias de dados da Residência 1'.format(app))
In [20]:
from keras.layers.core import Dense, Activation, Dropout
from keras.layers.recurrent import LSTM
from keras.models import Sequential
from keras.callbacks import ModelCheckpoint
from keras.models import load_model
from keras.optimizers import Adam
from keras.regularizers import l2
from keras.utils import plot_model
In [21]:
def build_fc_model(layers):
fc_model = Sequential()
for i in range(len(layers)-1):
fc_model.add( Dense(input_dim=layers[i], output_dim= layers[i+1]) )#, W_regularizer=l2(0.1)) )
fc_model.add( Dropout(0.5) )
if i < (len(layers) - 2):
fc_model.add( Activation('relu') )
fc_model.summary()
plot_model(fc_model)
return fc_model
fc_model_1 = build_fc_model([2, 256, 512, 1024, 1])
In [22]:
adam = Adam(lr = 1e-5)
fc_model_1.compile(loss='mean_squared_error', optimizer=adam)
start = time.time()
model_path = "./resources/ann-fc_refrig_h1_2.hdf5"
checkpointer = ModelCheckpoint(filepath="./resources/ann-fc_refrig_h1_2.hdf5", verbose=0, save_best_only=True)
hist_fc_1 = fc_model_1.fit( X_train1, y_train1,
batch_size=512, verbose=1, nb_epoch=200,
validation_split=0.33, callbacks=[checkpointer])
print('Tempo total de treinamento do modelo (s):', round(time.time() - start, 0))
In [24]:
fc_model_1 = load_model(model_path)
pred_fc_1 = fc_model_1.predict(X_test1).reshape(-1)
mse_loss_fc_1 = mse_loss(pred_fc_1, y_test1)
mae_loss_fc_1 = mae_loss(pred_fc_1, y_test1)
print('MSE no conjunto de teste: ', mse_loss_fc_1)
print('MAE no conjunto de teste:', mae_loss_fc_1)
In [25]:
train_loss = hist_fc_1.history['loss']
val_loss = hist_fc_1.history['val_loss']
def plot_losses(train_loss, val_loss):
plt.rcParams["figure.figsize"] = [24,10]
plt.title('MSE dos conjuntos de treino e teste - Resid. 1')
plt.plot( range(len(train_loss)), train_loss, color = 'b', alpha = 0.6, label='loss (treino)' )
plt.plot( range(len( val_loss )), val_loss, color = 'r', alpha = 0.6, label='loss (validação)' )
plt.xlabel( 'época' )
plt.ylabel( 'loss' )
plt.legend()
plot_losses(train_loss, val_loss)
In [26]:
plot_each_app(df1_test, dates[1][17:], pred_fc_1, y_test1,
'Rede Neural FC: Real e Previsão nos 6 dias do Conjunto de Teste da Resid. 1', look_back = 50)
In [32]:
y_pred_fc_2 = fc_model_1.predict(X_2).reshape(-1)
mse_fc_2 = mse_loss(y_pred_fc_2, y_2)
mae_fc_2 = mae_loss(y_pred_fc_2, y_2)
print('MSE no conjunto de teste: ', mse_fc_2)
print('MAE no conjunto de teste: ', mae_fc_2)
In [33]:
plot_each_app(df[2], dates[2], y_pred_fc_2, y_2, 'Modelo RNA-FC para o Refrigerador: treinado na Resid. 1, prevendo na Resid. 2')
In [35]:
def process_data(df, dates, x_features, y_features, look_back = 50):
i = 0
for date in dates:
data = df.ix[date]
len_data = data.shape[0]
x = np.array([data[x_features].values[i:i+look_back]
for i in range(len_data - look_back) ]).reshape(-1,look_back, 2)
y = data[y_features].values[look_back:,:]
if i == 0:
X = x
Y = y
else:
X = np.append(X, x, axis=0)
Y = np.append(Y, y, axis=0)
i += 1
return X,Y
In [36]:
start = time.time()
X_train, y_train = process_data(df[1], dates[1][:17], ['mains_1','mains_2'], df[1].columns.values[2:])
X_test, y_test = process_data(df[1], dates[1][17:], ['mains_1','mains_2'], df[1].columns.values[2:])
print('Tempo de execução total (s): ', time.time() - start)
print(X_train.shape, y_train.shape, X_test.shape, y_test.shape)
In [37]:
fc_model = build_fc_model([100, 256, 512, 1024, 1])
In [38]:
# Let's flatten data to feed into fc model
X_train_fc = X_train.reshape(-1, 100)
y_train_fc = y_train[:,2]
print(X_train_fc.shape, y_train_fc.shape)
In [39]:
adam = Adam(lr = 1e-5)
fc_model.compile(loss='mean_squared_error', optimizer=adam)
start = time.time()
checkpointer = ModelCheckpoint(filepath=model_path, verbose=0, save_best_only=True)
hist_fc2 = fc_model.fit( X_train_fc, y_train_fc,
batch_size=512, verbose=1, nb_epoch= 200,
validation_split=0.33, callbacks=[checkpointer])
print('Tempo total de execução (retreino - em segundos): ', time.time() - start)
In [40]:
train_loss = hist_fc2.history['loss']
val_loss = hist_fc2.history['val_loss']
plot_losses(train_loss, val_loss)
In [41]:
fc_model = load_model(model_path)
pred_fc = fc_model.predict(X_test.reshape(-1, 100)).reshape(-1)
print(pred_fc.shape)
In [42]:
mse_loss_fc = mse_loss(pred_fc, y_test[:,2])
mae_loss_fc = mae_loss(pred_fc, y_test[:,2])
print('MSE no conjunto de teste: ', mse_loss_fc)
print('MAE no conjunto de teste: ', mae_loss_fc)
In [43]:
plot_each_app(df1_test, dates[1][17:], pred_fc, y_test[:,2],
'Real e Previsto para o Refrigerator nos 6 dias do conjunto de testes da Resid. 1', look_back = 50)
In [44]:
start = time.time()
X_2, y_2 = process_data(df[2], dates[2], ['mains_2','mains_1'], ['refrigerator_9'])
X_2_fc = X_2.reshape(-1, 100)
y_2 = y_2.reshape(-1)
print('Tempo total de execução (s): ', time.time() - start)
print(X_2_fc.shape, y_2.shape)
In [45]:
pred_fc_50_h2 = fc_model.predict(X_2_fc).reshape(-1)
mse_loss_fc_50_2 = mse_loss(pred_fc_50_h2, y_2)
mae_loss_fc_50_2 = mae_loss(pred_fc_50_h2, y_2)
print('MSE no conjunto de teste: ', mse_loss_fc_50_2)
print('MAE no conjunto de teste: ', mae_loss_fc_50_2)
In [46]:
plot_each_app(df[2], dates[2], pred_fc_50_h2, y_2, 'RNA-FC retreinada (50 últimos dias): treinada na Resid. 1, prevendo na Resid. 2')
In [47]:
def build_lstm_model(layers):
model = Sequential()
for i in range(len(layers) - 2):
model.add(LSTM(
input_dim=layers[i],
output_dim=layers[i+1],
#stateful=True,
return_sequences= True if i < len(layers) - 3 else False ))
model.add(Dropout(0.3))
model.add(Dense(layers[-1]))
model.summary()
plot_model(model)
return model
model = build_lstm_model([2,64,128,256, 1])
In [54]:
# Habilitando HPU
import tensorflow as tf
if tf.test.is_gpu_available(cuda_only=False, min_cuda_compute_capability=None):
print('Ativando GPU!')
config = tf.ConfigProto(device_count = {'GPU': 1})
sess = tf.Session(config=config)
Out[54]:
In [56]:
start = time.time()
adam = Adam(lr = 5e-5)
lstm_model_path = "./resources/lstm_model.hdf5"
model.compile(loss='mean_squared_error', optimizer=adam)
checkpointer = ModelCheckpoint(filepath=lstm_model_path, verbose=0, save_best_only=True)
hist_lstm = model.fit(
X_train,
y_train[:,2],
batch_size=512,
verbose=1,
nb_epoch=200,
validation_split=0.3,
callbacks=[checkpointer])
print('Tempo de treino (s): ', time.time() - start)
In [57]:
train_loss = hist_lstm.history['loss']
val_loss = hist_lstm.history['val_loss']
plot_losses(train_loss, val_loss)
In [58]:
model = load_model(lstm_model_path)
pred_lstm = model.predict(X_test).reshape(-1)
print(pred_lstm.shape)
In [59]:
mse_loss_lstm = mse_loss(pred_lstm, y_test[:,2])
mae_loss_lstm = mae_loss(pred_lstm, y_test[:,2])
print('MSE no conjunto de teste: ', mse_loss_lstm)
print('MAE no conjunto de teste: ', mae_loss_lstm)
In [60]:
plot_each_app(df1_test, dates[1][17:], pred_lstm, y_test[:,2],
'Real e Previsto para os 6 dias do conjunto de testes da Resid. 1', look_back = 50)
In [61]:
pred_lstm_h2 = model.predict(X_2).reshape(-1)
mse_loss_lstm_h2 = mse_loss(pred_lstm_h2, y_2)
mae_loss_lstm_h2 = mae_loss(pred_lstm_h2, y_2)
print('MSE no conjunto de teste: ', mse_loss_lstm_h2)
print('MAE no conjunto de teste: ', mae_loss_lstm_h2)
In [62]:
plot_each_app(df[2], dates[2], pred_lstm_h2, y_2,
'Real e Previsto para os 6 dias do conjunto de testes da Resid. 2', look_back = 50)
In [ ]:
In [ ]: